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ABSTRACT 



Finite Markov processes are considered, with bi -dimensional state 
space, such that transitions from state (n,i) to state (m,j) are possible 
only if m < n+1. The analysis leads to efficient computational algorithms, 
to determine the stationary probability distribution, and moments of first 
passage times. 



1. INTRODUCTION 



We consider irreducible, bi-dimensional Markov processes 
{(X(t) ,Y(t)) , t > 0), on a finite state space {(n,i), 0 < n < N, 1 < i < K^}, 
for which transitions from state (n,i) to state (m,j) are possible only if 
m < n+1. The infinitesimal generator for such a process has a lower block- 
Hessenberg form : 



o 

o 

O' 


^ 0,1 


0 


0 


... 0 


0 


^ 1,0 


‘^l.l 


<^ 1.2 


0 


0 


0 


^2,0 


^ 2,1 


^2,2 


^2,3 


... 0 


0 


^N-1,0 


^N-1,1 


Vl ,2 


^N-1,3 


”• ^N-1,N-1 


^N-1,N 


^N,N 


^N,l 


^N,2 


'^N,3 


••• *^N,N-1 


'^N,N 



The diagonal blocks ^ are square matrices of order K^, 0 < n < N. 
Their diagonal elements are strictly negative, all other elements are non- 
negative. The blocks Q , n + m, are rectangular with appropriate dimen- 

n jfii 

sions; their entries are non-negative. The rowsums of Q are equal to zero. 

We denote by level n the set {(n,i), 1 < i < K^} of states corres- 
ponding to the common value n for the first index. The structure (1.1) of 
Q permits the Markov process to move up by only one level at a time, and to 
move down by several levels at a time. The process is said to be skip-free 
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Many models in queueing theory lead to Markov processes (or Markov 
chains) for which the infinitesimal generator (respectively the transition 
probability matrix) is a (possibly infinite) block-Hessenberg matrix. 

The paradigm for the lower-Hessenberg case is the GI/M/1 queue, for the 
upper-Hessenberg case, it is the M/G/1 queue. 

Neuts [7,9] has studied in great detail Markov processes and Markov 
chains with infinite block-Hessenberg matrices which exhibit a repetitive 
structure : for some n*, 

'‘n = '‘r. • 

for all n>n», n,<nn, 

if the matrix is lower-Hessenberg, and a similar condition if the matrix is 
upper-Hessenberg . 

In particular, Neuts has shown that if a process with such a lower 
block-Hessenberg matrix is positive recurrent, its stationary probability 
vector has a modified matrix-geometric form. For processes with upper 
block-Hessenberg matrix, the quantities which may be most elegantly studied 
are first passage times to lower levels; these first passage times correspond 
to the busy period in queueing theory. 

There do not exist many published results for processes with f i ni te 
block-Hessenberg matrices. Most of them deal with block-Jacobi matrices (the 
blocks are zero if either m > n+1 or m < n-1). Torrez [12] proposes to 
determine the stationary distribution for birth-and-death processes in a random 
environment, by using numerical procedures designed for band matrices. 

Neuts [10] and Carrol et al . [1] explicitely determine the stationary probability 
distribution for the PH/M/C and M/PH/1 queues with finite waiting room. 

Hajek [3] considers a model with repetitive structure •’ Q„ _ = Qi ,, 1 < n < N-1, 
Qn "'^12 ^n+1 n " ^2 1’ 1 ^ ^ ^ ^"2- Processes with general block- 

Jacobi matrices are analysed in Keilson et al . [4], and in Gaver et al . [2], 
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where specific numerical examples are also given. 

Wikarski [13] determines, by purely analytic arguments, the stationary 
probability vector for processes such that Kq ^ ^ ••• ^ such that 

the matrices are of rank 0 < n < N-1. We do not need these 

restrictions. Moreover, our approach yields quantities which have a clear 
probabilistic interpretation. 

The purpose of the present paper is to show how the stationary probabi- 
lity distribution, and moments of first passage times from some level to another 
level, may be efficiently computed. 

The stationary probability distribution is analysed in the next section. 
In Section 3, we analyse the first passage time from level n to a higher level 
m, 0 < n < m < N. The results in these two sections are related to those in [2]. 

The analysis of the first passage time from level n to a lower level k, 
0<k < n<N, is more involved, and is presented in Section 4. We shall 
observe that it is possible to compute independently the stationary probability 
distribution, and the moments of first passage times to higher levels. On the 
other hand, it is necessary to compute both the stationary distribution, and 
moments of first passage times to higher levels, in order to compute the mo- 
ments of first passage times to lower levels. 

In Section 6, we give some numerical examples of common-cause of 
failure models. We consider systems with N machines and R repairmen. If n 
denotes the number of machines in working condition, the repair rate is 
umin(R, N-n). The failures occur in this way : in an interval (t, t+dt), there 
is a jump from n to n-1 with probability xn dt + o(dt) (independent failures) 
or to n-k, 1 <k < n, with probability v([^)q*^( l-q)’’"*^dt + o(dt) (common-cause 
failures) . 

The results to follow are applicable, mutatis mutandis, to Markov pro- 
cesses with upper block-Hessenberg infinitesimal generator, and to finite 
Markov chains with block-Hessenberg transition probability matrix, as we briefly 



mention in Section 5. 
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Notational conventions . 

We define a number of vectors and matrices with indices. Occasionally, 

we refer to specific elements of those. In order to somewhat simplify the 

notations, we use Q_ „i(i,j) to represent the (i,j)th element of the matrix 
n j n 

*^n,n' ■ 

We systematically use "n" to denote a given level, "m" to denote a 
higher level, and “k" to denote a lower level. 

2. THE STATIONARY PROBABILITY DISTRIBUTION 



We denote by £ the vector of stationary probabilities associated with 
the Markov process Q : 



£Q=0, £6=1, 



( 2 . 1 ) 



where £ is a column vector with unit elements. We partition the vector £ as 
£ = (Pq, £j^, ..., £j^), where the subvectors have elements and correspond 
to the states of level n, 0 < n < N. 

In order to determine the vector £, we need three preliminary lemmas. 

Lemma 1 

Consider a Markov process on the state space {1,2, . . . ,r,r+l,r+2, . . . ,r+s}, 
with infinitesimal generator Q : 



'V 

Q = 



A A 

0 0 



where A is a square matrix of order r, A is a rectangular r by s matrix. 
The states r+1 to r+s are all absorbing. 

a) The states 1 to r are all transient if and only if the matrix A is non- 
singular, in which case (-A"^) > 0. 
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b) The (i,j)th entry of (-A~^), for 1 < i , j < r, is the expected amount 

of time spent in the transient state j, starting from the transient state i, 
before absorption in any of the absorbing states. 

c) The (i,k)th entry of (-A"^ A), for 1 < i < r, r+1 < k < r+s, is the probabi- 
lity that, starting from the transient state i, absorption occurs in the 
state k. 

d) The (i ,k)th entry of (-A'^)(-A"^ A) , for 1 < i < r, r+1 < k < r+s, is the 
expectation of the time spent in the transient states, starting from the 
transient state i, before absorption, measured on those paths that lead to 
absorption in the state k. 

Proof . These results are well known, although we are not aware of publications 
where they are proved under the stated form. The first statement is proved in 
Neuts [9], Lemma 2.2.1; the second statement is a consequence of Neuts and 
Meier [11], Corollary 2; the third statement is proved in [2], Lemma 1. 

To prove the fourth statement, we must show that the (i,k)th entry of (-A"^)* 
("A"^ a) is equal to t d v^ k^^^’ 

'^■i i-(t) ■ P process gets absorbed in k, before time t 1 the 
initial state is i]. 

Conditioning on the state of the process at time t > 0, one obtains, for h > 0 
Vj_|^(Wh) = (1*A„ h)v,_|^(t) * h * Z A,J h * 0 (h), 

J ▼ I 

hence vl^^(t) = Z A.j Vj^^(t) + A.^^ , 

J 

or v'(t) = A v(t) + A , 

where v(t) is the matrix formed by the v,j j,(t) 's. We already know that 
v(0) = 0, and lim v(t) = (-A'^)A, therefore, we easily obtain that 
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/” t d v(t) = - [t((-A"^)A - v(t))]J 



+ ({-A‘^)A - v(t))dt 

= Q {-A’^) V(t) dt 
= (-A‘^) [v{t)]J = (-A)'^-A'^A). 



□ 



We now define the processes S^, 0 < n < N. For 1 < n < N, is 
the restriction of the original process Q, observed during those intervals 
of time spent at level n, before the process Q enters level n-1 for the 
first time. The state space of is {(n,j), 0 < j Clearly, all 

§n. 1 ‘^n <N, are transient Markov processes. The process Sq is the res- 
triction of the process Q observed at the lowest level 0; it is an ergodic 
Markov process. We denote by the infinitesimal generator of the process 
S^, 0 < n < N. 

We also define the matrices n ,,0<k<n<Nas follows : 

n J K 

is the probability that, starting from state (n,i), the first 
n j K 

state visited at any level below level n is the state (k,j). 

Lemma 2 

The matrices C^, 0 < n < N are recursively determined as follows ; 



Proof . The equation (2.2) is obvious : starting from level N, the process 
S|^ terminates as soon as the process Q enters any of the levels 0 to N-1. 
Meanwhile, the transitions are governed by 

Consider n fixed, 1 <n <N-1, and let Z^(t), t > 0 be defined as follows. 
Zn(T) = i if the process is in state (n,i) at time x. 



" ^N,N ’ 



( 2 . 2 ) 




, 0 < n < N-1. 



(2.3) 
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Assume that Z^{x) = i and Z^(t + dt) = j + i. Since j + i, a transition 
must have occured in the process Q. Z^(t + dx) = j if and only if one of 
the following events have occured. 

a. The transition is from (n,i) to (n,j), this happens with probability 
Qn,n dT + o(dx). 

b. The transition is from (n,i) to (n+l,h) for some h, and the process Q 

returns to (n,j), after spending an unspecified amount of time at level 
n+1 or higher, before visiting any other state at level n or lower. This 
happens with probability Qn^n+x(i.h) dx + o(dx). 

It is now clear that 



which proves (2.3) for n > 1. The proof for n = 0 is identical. 

Lemma 3 

The matrices n ,,0<k<n<N are recursively determined as 
n y K 



P[Z^{x + dx) = j1Z„(x) = i] 




+ o(dx), for 1 < i * j < K^, 



and similarly, that 



P[Z^(x + dx) = i|Z„{x) = i] 




+ o(dx) , for 1< i < K^, 




for 0 < k < n < N-1 (2.5) 



(2.4) 



fol 1 ows 
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Proof . The equation (2.4) is an immediate consequence of Lemma 1 and (2.2). 
A simple probabilistic argument leads to 



Remark . The matrices C^, 0 < n < N and n^|^,0<k<n<N, must 
be numerically computed together, as indicated below. 

Algorithm A 



A. 2 for k from 0 ^ N-1 ^ determine rTj^ from (2.4) ; 



begi n determine from (2.3); 

for k from 0 to n-1 do determine n_ ^ from (2.5) 

' ' ' ■ - - — - n > K 

end ; 

A. 4 determine Cq from (2.3)* 

We are now in a position to prove our first main result. 
Theorem 1 . 

The vectors 0 < n < N, are determined by the equations 




which, together with (2.3), implies (2.5). 



□ 




A. 3 



fo r n from N-1 1 step -1 do 



= Vl Vl,n for 1 <n<N, 



( 2 . 6 ) 

(2.7) 



N 




( 2 . 8 ) 
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Proof . From the structure (1.1) of Q, it results that the system £ Q = 0 
may be decomposed into 

£o<lo.o‘ ", Ev<lv. 0 -i- (2-3) 

V=1 

£k-l «k-l,k * Sk \.k * K «..k • “• f"- 1 <"-l. (2-W) 

En- 1 '^N-l.N * En *^«,N ° (2-11) 

From (2.11) and Lemma 2. £^ = -£„.^ (-Q-J^). 

which proves (2.7) for n = N. 

The remainder of the proof uses an induction argument in two steps. 
Induction assumption 1 : for a given k, 1 < k < N-1, the equation (2.7) 
holds for k+1 < n < N. 

In order to prove that (2.7) holds for n = k, we first need to show 

that 

«.,k = "t.k • f'”' (2->2) 

Clearly, (2.12) holds for t = N, by (2.4). An induction argument will be 
used to prove (2.12) for arbitrary t > k+1. 

Induction assumption 2 ; for a given a, k+1 < o+l < N, the equation (2.12) 
holds for a+1 < t < N. 

We prove that it holds for t = o, by using successively the assump- 
tion 2, the assumption 1, and (2.5). 

N 

^ Hv ^v.k " fia ^a.k £a+l ^"^a+l^ Vl.k ’ 
v=a 

* E« (‘^a.k * '*a,a-U "a*l,k)’ 

■ E» <-'«) "a.k • 
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Finally, to show that (2.7) is satisfied for n = k, note that 
(2.10) can be written as 



N 

- = Hk-1 ^k-l,k Hk ^k.k £v ^v.k 

' £k-l ^k-l,k * fik “^k.k ■*■ £k+l ("^k+l^ ^k+l,k 
"" Hk-i Qk-l,k ■*■ £k (^k.k ■*■ ^k,k+l '^k+l.k^' 



The equation (2.6) is proved in the same way. Since Cq is the 
infinitesimal generator of a finite, irreducible Markov process, the equation 
(2.6) has a unique solution, up to a multiplicative constant, and that 
constant is determined by (2.8). a 

We have thus shown that the matrices C^, 0 < n < N, are the 
crucial quantities needed to determine the stationary probability distribution. 
We shall observe in Section 4 that the same matrices, together with the 
matrices 0 < k < n < N, play a central role in the evaluation of first 

passage times from a given level n to a lower level k. We must, however, 
first examine the first passage times to higher levels. 



3. FIRST PASSAGE TIME TO HIGHER LEVELS 

We denote by t • the first passage time from level n to a different 

I I 9 I I 

level n' : = inf {t > 0 : X(t) = n'|X(0) = n}, 0 < n,n' <N, n + n*. 

We are concerned in this section with the first passage times t for 

n 

m > n. We need to first establish some preliminary results 

We define the processes S^, 0 < n < N in a similar fashion to the 
processes S^. For 0 <n < N-1, the transient Markov process is the restric- 
tion of the process Q, observed during those intervals of time spent at level n, 
before the process Q moves upward to level n+1 for the first time. The ergodic 
process is the restriction of the process Q observed at the highest level N. 
We denote by the infinitesimal generator of the process S^, 0 <n <N. 
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We also define the matrices 0^n<m<Nas follows : 

the probability that, starting from state (n,i), the first state visited at 
level m is the state (m,j). More formally, r„ „(i.j) = ) = jlX(O) = 

n, Y(0) = i]. 

Clearly, I’n ^ — “ — ’ (3-1) 



Lemma 4. 



The matrices C^, 0 < n < N, are recursively determined as follows 



^0 ° '^ 0,0 ’ 



The proof is analogous to that of Lemma 2 and is omitted here. 



(3.2) 

(3.3) 



Lemma 5 

The matrices r „,0<n<m<N are recursively determined as 
n,m 

fol 1 ows 

^n,m ^n,m-l ^m“l,m (3.5) 

m 

= n r. , , for 0 < n < N-2, (3.6) 

t=n+l 

n+2 < m < N . 



Proof . For n=0, the equation (3.4) is an immediate consequence of Lemma 1 
and (3.2), For 1 < n < N-1, observing the system at the first departure from 
level n leads to 

-1 

^n,n+l ^'^n,n^^*^n,n+l ^ ^n,k ^k,n ^n,n+l^’ 

which, together with (3.3), implies (3.4). To prove (3.5), since Q has the 

structure (1.1), one nnerely decompose the probabilities in r according to 

n >rn 
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the first state visited at level m-1. Equation (3.6) is now obvious. 

□ 

The matrices C„, 0 < n < N and r„„, 0<n<m<N, must be numeri- 
n n,m 

cally computed together, as indicated below. 

Algorithm B 



B.l Cg := Qq^q ; 

determine Fg ^ from (3.4); 

B.2 fo^r m from 1 jto N-1 do 

begin determine from (3.3); 

determi ne from (3.4); 

for n from 0 m-1 ^ determi ne from (3.5) 

end; 

B.3 determine Cj^ from (3.3). 



We now define, for x>0, 0<n,n' <N, n* n', l<i 
1 < j <K^. , 



9n,i ;n' 



,j(x) = <x. 



Y(Tn,n-+ 



we denote by G „i(C) the matrix such that G _i(C)(i.j) equals the Laplace- 
n j n n j n 

Stieltjes transform of the possibly defective distribution g . • •(.). 

n > 1 *n I j 

Clearly, since the Markov process may only move up by one level at 
a time, we have that 



m 

= n G. , .(C), for 0 < n < N-2, 
t=n+l 

n+2 < m < N. 



(3.8) 



(3.9) 
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A simple argument conditioning on the state of the process at time h > 0 
shows that 






v+l 
n-1 



* k=0 ji*l 

^ <Qn,n+l)iJ forl<n<N-l. 



(X) 



where ♦ denotes the Stieltjes convolution, from which we conclude that 



^ *^n,n+l^^^ “ *^n,n ®n,n+l^^^ ^n,n+l 



n-1 



G„,„,l(e).forl<n<N-l. 



and similarly, 

^ S.l(^) = ^0.0 ^ %.r 

This equations may be written as 

Gn,n.l(^) = ^n,n.l * for 0 < n < N-1, 

where Dq(?) = (cl - Qq.o)'^ * 

»r(E) = - «r.n ' '’"■k '«'• * <"-1. 

Of course, we have jjj(O) = » sod from (3.3) D^(0) = (-C”^), 

0 < n < m < N. 



(3.10) 

(3.11) 

(3.12) 

(3.13) 

(3.14) 
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d 



Let 


K n* = 
n ,n 


■ 3 ^ k=0 * 




(3.15) 


and 


= 

— n,n 


'^n.n' - ’ for 0 < n,n' 


< N, n * n' . 


(3.16) 


We have 




i) = E[t„ Y(t„ „ + 0) = 

" n,m ^ n,m ' 


j|X(0) = n, Y(0) = i], 




and 




= I X(0) = n. Y(0) 


= i]. 





Theorem 2 . 

The expected values of the first passage times to higher levels 
are given by the following relations. 



^,1 " ^ 



0 > -» 






n-1 



'in.k 1< n < N-1. 



y„™=u„„,+r„_,u_i_, for 0 < n < N-2, 
•m»in — n^ni”! n>ni“l — m-ijni 

n+2 < m < N . 



(3.17) 

(3.18) 

(3.19) 



Proof . The equation (3.17) immediately results from Lemma 1. The 
equations (3.19) are directly proved by using the strong Markov property : 

^ I m - V(0) . 11 

° ^ I " 

Conditioning on the first level visited after level n, and using Lemma 1, 



we obtain 
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'n.k ^,n+l ’ 




n,k ^-^.n ^ ^k,n ^^,n+l^’ 



by (3.19). Hence, 



n-1 



n-1 




which, together with (3.3), proves (3.18). 



□ 



Remarks. 



a. Another argument, using the equations (3.12), (3.15) and (3.16), 
requires the evaluation of the matrix 3 Dj^(c)/3C . This is done in 
Appendix A. Similarly, higher moments of first passage times may be 
obtained by purely routine, albeit tedious, manipulations. 

b. The vectors u^ may be computed together with the matrices by 

appropriately modifying Algorithm B. 

4. FIRST PASSAGE TIME TO LOWER LEVELS 

We denote by the first passage time from level n to a lower 



level k. In the preceding section, we have in essence decomposed the first 

passage time from level n to a higher level m, into a sum of first passage 

times from n to n+1, from n+1 to n+2,..., from m-1 to m. Since it is 

possible to move down several levels in one transition, no such simple 

decomposition is possible for x . . Rather, we must explicitely consider 

n ^ K 

the first level visited below n, starting from level n. In order to do this, 
we define the first passage time T^ as follows 



Tn = inf {t > 0 ; X(t) < n | X(0) = n} , 

= min {x k . 0 < k < n-1}, for 1 < n < W. 
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We define, for x>0,0<k<n<N, l<i<K^, 

9n,i;k.j('') = P + 0) = k. Y(T^ + 0) = jlX(O) = n, Y(0) = i], 

we denote by G i,(C)(i>j) the Laplace-Stiel tjes transform of the possibly 
n 9 K 

defective distribution g ... •(.), and by G .(C) the matrix with entries 

itylyKyJ N>K 

given by Gj^^k(5)(i ,j). 

An argument similar to the one used to obtain (3.10) and (3.11) shows 

that 



^ ‘^N.k * *^N,N for 0 < k < N-1 



(4.1) 



and c = Qn.k + ^n.n ^n.k^^) 



(4.2) 



^ ^n.n+1 ^^n+l.k^^) ^n..i n(0 K for 0 < k < n < fl-1. 



n+l,n'^' n,k 



These equations may be written as 



GN,k(0 = ^N,k ’ forO<k<N-l, 



(4.3) 



Sn.k(^) = [Qn.k " Qn,n.l ^n.l.k^^)^* for 0 < k < n < N-1 . (4.4) 



where Djj(c) = (cl - Qn w)"^ , 



N.N^ 



-1 



^n(^) = - Qn.n ' ^n.n.l S.l.n^^)) * for 1< n < N-1. 



(4.5) 

(4.6) 



Obviously, we have that 6_ . (0) = n . , and (Lemma 2) D (0) = (-C"^) 

n>K iijK n n 

Observe that 



n . e <e, for 0 < k < n < N, 

n ) K 

n-1 

Z n . e = e, for 1 < n < N. 



k=0 



n,k — 



(4.7) 
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^n,k -- ^5=0 



(4.8) 



3C 



On^k('*J) = = I'* + 0) = j|X(0) = n. Y(0) = i], 



Lemma 6. 



The matrices 0^^ are recursively determined as follows 
^N,k " ’^N.k 



k = (-Cn^) n. . , for 0 <k <N-1, 



(4.9) 



^n.k ■ ("^n ■*■ n+l(^n+1 k ■*■ ^n+1 n ^n.k)^* 



^n.k - ^‘"n.k ^ ''n,n+l^"n+l,k ^ "n+l,n "n.k^'* 

for 0 < k < n < N-1. 

Proof . It results from Lemma 1 that 
'"1\ / ^"1 



^N,k " ^N,k^’ forO<k<N-l, 



which together with (2.4) proves (4.9). For n <N-1, we condition on 
the first level visited after level n, and we obtain by using Lemma 1 and 
the strong Markov property that 

^n,k “ ^"*^n,n^ *^n,k ^~^n,n *^n,n+l^^^n+l,k ^ ^n+l,n *^n,k *^n+l,n ^n,k^* 



hence 



^*^n,n *^n,n+l *^n+l,n^ ^n,k ^*^n,k ^n,n+l ^^n+l,k ^n+l,n 



which, together with (2.3), proves (4.10), 



Corollary 1 . 



The vectors u„ defined as 
— n 

n-1 _ 

^ = Z 0„ . e 
-n n,k - 



(4.11) 



are recursively determined as follows 



■ (-c;;') i . (4.12) 

. forl<n<N-l. (4.13) 

□ 

This corollary immediately results from Lemma 6 and (4.7). 

The Laplace-Stiel tjes transforms of the distribution functions 
for the first passage times t ^ are given by the matrices G ^(C). If 
one decomposes the trajectories from level n to level k, according to the 
first level v visited below n, one immediately obtains 



k-1 



“"n.k^^) ,.(5) G,.k(C) 



n,k 



v=0 



n,v' 






n-1 

\xi (C) G .(c), forO<k<n<N, 

v=k+l n,v' ' v,k'^' ’ 



(4.14) 



from which the moments of first passage times to specific lower levels 
may be derived. The first moments may be more easily obtained by using 
the strong Markov property. The next theorem is stated without proof. 



Theorem 3. 



The expected values of the first passage times to lower levels 
are given by the following relations 



^1.0 


= Hr 










k-1 


n-1 




= + 


z n 


u , + Z 


-HTI y k 


— n 


v=0 


v=k+l 



(4.15) 

for 0<k<n<N, (4.16) 



where the matrices n are determined by Lemma 3, and the vectors u . 

n > V ~“V ) K 

for V < k-1 are determined by theorem 2. □ 
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We observe here a major difference between the first passage times 
to higher and to lower levels. We have observed in the preceding section 
(Remark b. after Theorem 2) that the vectors for m > n may be computed 
in one loop, starting from n = 0, to n = N-1. The evaluation of the vectors 
1 ^ 1 ^ for k < n requires two steps. One for computing the vectors starting 
from n = N, to n = 1. In the second step are computed the vectors u^ 
starting from n = 1, to n = N. We observe also that the vectors ^ for 
m > n have to be computed, even if one is only interested in the first passage 
times to lower levels. 

It is clear that one can combine the computation of the stationary 
probability distribution, and the expected first passage times, as is done in 
the algorithm C given in Appendix B. 

It is much easier to obtain moments for the first passage time from 
level n to any level at or below k < n-1. Let 

Tn,k = inf {t > 0 : X{t) <k|X(0) = n}, 

= min {t : 0 < V < k} , for 0 < k < n < N. 
n j V 

We denote by G .(c) the matrices with entries equal to the Laplace-Stieltjes 
n j K 

transforms of the distribution functions for the first passage times T^ 

Also, let 

and I, = U„ 1 , e , for 0 < k < n < N. 

— n,k n,k — 



If we condition on the first level v visited below n, we obtain 



n-1 



h,k 



(C) = I G (c) + I G (C) G |^(C). 



v=k+l 



from which the moments of T^ may be obtained. Again, it is easier to deter- 
mine the first moments by using the strong Markov property. The next theorem 
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is stated without proof. 



Theorem 4. 



The expected values of the first passage time from level n to any 
level at or below k, are given by 



and 



^,n-l ■ -Hn ’ 



^ n 

u„ . = u„ + z n„ u . , 
— n,k n y=k+i — v»k 



for 0 < k < n < N, 



(4.17) 

(4.18) 



where ^ is defined in (4.11). □ 

The vectors ^ may be computed in one loop, beginning with 
k = N-1, n = N, and progressively decreasing k until k = 0; for each value 
of k, one computes the vectors ^ to |^. In Section 5, we shall refer 
to Algorithm 0 (in Appendix B), which is designed to recursively compute the 
vectors |^, N-1 >k >0. 



5. GENERAL REMARKS 



Remark a. 



We firstly comment on Equation (2.7). The (i,j)th entry of the 
matrix ^ equal to 

‘Vl.n 
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The factor probability that, 

upon leaving the state (n-l,i), the process Q moves to the state (n,k). 

The (k,j)th entry of is the expected time spent by the process in 

state (n,j), starting from (n,k), before hitting any state at a lower level, 
by Lemma 1 and the definition of If that lower level is below level n-1, 
the process will be forced, by the structure (1.1) of Q, to pass through level 
n-1 before reaching again level n. 

Therefore, it results from (5.1) that [Q , „(-C"^)](i ,j) Is equal to 
(-Qn_i p_j)(i »i)times the expected time spent in the state (n,j), before the 
next return to level n-1, given that the process Q starts in the state (n-l,i). 
This is exactly the interpretation of Neuts' matrix R ([9], Theorem 1.7.2, and 
[8]), and shows that (2.7) is the counterpart to the matrix-geometric solutions. 

Remark b . 

If the matrix Q is upper block-Hessenberg, instead of having the form 
(1.1), then our results can be modified in an obvious way. In that case, the 
first passage times to lower levels may be directly computed, via an algorithm 
derived from Algorithm B, while the first passage times to higher levels are 
more involved. Also, the stationary probability vector is computed by determi- 
ning first the vector (up to a multiplicative constant) and then recursively 
the vectors .... £q. 

This last observation points to the reason why it is difficult to eva- 
luate the stationary probability distribution for infinite upper block-Hessenberg 
matrices. 

We briefly examine now the problem of obtaining a good approximation for 
the stationary distribution in such a case. Neuts [7] and Lucantoni and Neuts 
[6] consider infinite upper block-Hessenberg matrices of the following form ; 
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B, 



0 



B 



1 



B, 



’2 



b 



3 



B, 



0 



A 



1 



A. 



'2 



A 



3 




0 0 Ajj Aj 



and stationary probability vector x = ()^, x^, An explicit form 

is found in [6,7] for and Xj^, as well as complicated but tractable expres- 
sions for 



From these, a heuristic method is proposed to truncate the matrix Q at some 
level N, and it is suggested that x^, ^ be determined by a block 

Gauss-Seidel iterative procedure. 



1, the vector £ corresponding to the finite matrix, and to compare the vectors 
{£o»£_i) to in order to check whether the level N is properly chosen. 



Also, we may truncate the matrix Q on the basis of first passage times 



argument, instead of using m^ and m 2 . It is very easy to adapt Algorithm D 
(in Appendix B), in order to compute successively the expected first passage 

time from level 0 to any level at or above m, for m = 1,2 One may then 

choose N to be the smallest value of m for which the expected first passage 
time exceeds a given, large value. When N is determined, the vector £ can be 
readily obtained, since most of the preliminary computations are already 
performed. 



oo 



m, = Z V X e 
1 ^ - 




v=0 



An alternative method would be to determine, as indicated in Section 
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Remark c . 

We have not considered first passage times to specific states, but 
have concentrated our attention on first passage times to levels. The reason 
is that, in a number of applications, such as chains in a random environment, 
the levels themselves are of importance. Nevertheless, it is not very diffi- 
cult to extend our analysis to the study of t ^ .„ ,• = inf{t > 0 : X(t) = m, 
Y{t) = j I X(0) = n, Y(0) = i}. With the results of Sections 3 and 4, we can 
observe the process between successive visits to level m, until, for the first 
time, the state (m,j) is visited before the process leaves level m. 

Remark d . 

We have only considered continuous-parameter Markov processes. Our 
analysis may be immediately adapted to discrete-time Markov chains. The major 
difference is that, in all our equations, the matrices and will 

be replaced by fundamental matrices for absorbing Markov chains (Kemeny and 
Snell [5], Section 3.2). 

6. NUMERICAL EXAMPLES 



We have computed some characteristics for different machine-repairmen- 
like processes, using the concepts of common-cause failure (described in the 
introduction) and of randomly varying environments. We have considered four dif- 
ferent types of systems, described below in increasing order of their complexity. 

We have computed, for each system, the stationary distribution and the 
vectors ^ |^. To do this, we have used Algorithm E (in Appendix B), which inclu- 
des parts of both Algorithms C and D. 
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We shall not attempt to answer any specific question about such 
systems. Rather, we shall describe the qualitative and quantitative diffe- 
rences in their behaviour. We shall use here, as the main performance measure, 
the expected times u^ until there remain k or less operating machines, start- 
ing from a fully operational system. 

Type I . The first model is the classical machine- repairmen process. We con- 
sider a system of N machines and R repairmen, with failure and repair rates 
respectively given by x' and y. The infinitesimal generator is in this case 
a tri -diagonal matrix. 

Type II . This is a simple common-cause failure model. In superposition to 
a machine-repairmen system, with N machines, R repairmen, failure and repair 
rates given by x and y, we consider an independent point process with exponen- 
tial interevent time distribution, with parameter v : the process of common- 
cause failures. Whenever a common -cause failure occurs, each operating machine 
has a probability q of failing, independently of the other machines. Thus, if 
there are n operating machines, the number of machines that fail due to common- 
cause has a binomial distribution, with parameters n and q. The infinitesimal 
generator has the form (1.1), with blocks of order 1. 

We define as system I I. A the common-cause failure system with N = 10, 
R = 1, X = .005, y = 1, V = .01, q = .2. If 1/y is equal to one hour, then 
1/X «8 days, and 1/v «» 4 days. The system I. A is a machine-repairmen system 
with N = 10, R s 1, X' = .007, y = 1. That value for x’ is chosen so that both 
systems are globally similar in the following sense ; when n machines are ope- 
rating, they both have the same expected number of failures in an interval of 
length dt : this is n(X dt + q v dt) + o(dt) = .007 n dt for system II. A, and 
n x' dt + o(dt) = .007 n dt for system I. A. 
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We give in Table I the stationary probability distribution for these 
two systems. We observe that both systems have high, nearly equal, probabili- 
ties of having most machines in operation. Nevertheless, the two systems are 
different since the distribution II. A has a longer tail towards the states 
with few operating machines. 

This difference is strikingly apparent on Table II, where we give the 
expected times Ujq 0 < k < 9. For instance, if we compare Ujq ^ for both 
systems, we observe that the average time until 5 machines or more are under 

3 

repair, is smaller by a factor 10 for the common-cause failure system. 

We conclude that multiple failures have a profound effect on the dyna- 
mic behaviour of the system. In support of our conclusion that multiple failures 
are the main cause for difference, we observe that both systems have identical 
repair facilities (R = 1, w = 1), and that the stationary expected number of 
failures per unit of time (henceforth denote by cp) are nearly equal : they are 
respectively given by 

10 

‘Pt A * ^ Pn = .0694815, 

^ n=0 ” 



and 



10 

'^II A ” ^ ^ “ .0693535. 



Type III. We assume now that the failure and repair parameters change according 
to an independent environmental process : they are respectively given by X-, jj-, 

J J 

V. and q. when the environment state is j, 1 < j < K. We assume that the envi- 

J J 

ronmental process is an irreducible Markov process, with infinitesimal generator 
T. The system is now described by a pair of indices (n,j), where n denotes the 
number of operating machines, and j denotes the state of the environment. The 
infinitesimal generator has the form (1.1), with blocks of order K. A typical 
example, with N = 4, R = 2 and K = 3 is displayed in Appendix C. Note that the 
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off-diagonal blocks are diagonal matrices because the only possible transi- 
tions are from (n,j) to (n,j'), with j * j', or (n,j) to (n'.j), with n + n'. 

We firstly define a system with two environments. In environment 1 
(the off-environment), no common-cause failure may occur. In environment 2 
(the on-environment), common-cause failure may occur at the rate v*. Indivi- 
dual failures occur with the same rate X in both environments. 

In order to make meaningful comparisons with the system II. A, we have 
chosen v* so that the expected time between two conmion-cause failures is equal 
to 1/v = (.01)'^, the same as for the system II. A. The time between two succes- 
sive common-cause failure has a PH-distribution (Neuts [9], chapter 2), with 
representation (6 , V) , where 

[ ■^^ 21 ‘^ ^21 



6 = [1 0], V = 



and its first moment is given by -6 V 




12 



(T21+ T 12 )(Ti 2 v*)"^, hence 



V* = V { 1 + ) . 

^12 

We have chosen T^2 = .00015 and I 21 ~ -004, from which we obtain 
V* = .27666. If the unit of time is one hour, then 1/1^2'^ ^ months, I/T22 
10 days, l/v* « 3.5 hours. 

We now define the systems III. A and III.B. For system III. A, we have 
N = 10, R = 1, K = 2, Xj^ = X2 = .005, = 42 “ 1» = 0, V2 - .27666, q^^ = 0, 

q2 = .2. For system III.B, there is only one difference : U2 = 2, the service 
rate is doubled in the on-environment. 
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In Table III, for each system, the first column represents tne 
stationary marginal distribution of the number of operating machines; the 
last row represents the stationary marginal probability of being in each 
environment; in the remaining two columns are given the stationary conditional 
distributions of the number of operating machines, given that the system is in 
the corresponding environment. 

The comparison of the column III.A.O in Table III and the column II. A 
in Table I shows that the two distributions are very similar, but the marginal 
distribution for the system in a random environment has a longer tail towards 
the states with few operating machines. We also observe that the marginal 
distribution is not a good descriptor for the system in a random environment : 
the two conditional distributions are quite different, which may have important 
implications, even if there is a small probability of being in the on-environment. 

In Table IV are given, for the same two systems, the values of Uj^q |^(1) 
and Ujg 1^(2), for 0 < k < 9. The last row indicates the conditional probability 
of being in each environment, given that there are 10 machines in operation. 

The columns III.A.l and III. A. 2 are to be compared to the column II. A in 
Table II. 

Several differences may be observed. Firstly, the time u^g g to com- 

3 

plete failure is shorter by a factor 10 for the system in a random environment. 

Secondly, the values of Ujg ^(1) to Ujg ^(1) are close to each other ; they 

differ at most by a factor 2. This is due to the influence of the on-environment 

-1 3 

it takes on the average (T^2) ^ units of time to switch from environ- 

ment 1 to environment 2, at which time the common-cause failure process is 
turned on and brings more rapidly the system towards the states with few opera- 
ting machines. As Uj^g g(l) is so close to computed 

Ug 1^(2), for 0 < k < 5, and found that it differs from the corresponding Uj^g |^(1) 
by at most 70 units of time. 
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We conclude that the systems 1 1. A and 1 1 1. A have markedly different 
dynamic behaviours, even when they have identical repair facilities and 
nearly equal stationary failure rates : 

10 

'^III.A " " ^Pn,l ^1 * Pn,2 ^^2 * “ .0666435, 

is close to <Pj| 

When examining the data for system III.B, we observe that the same 

comparisons can be made between III.B and II. A as between III. A and II. A. 

10 

The stationary average service rate for system III.B is E (p„ i + p„ o mo) = 

n=0 ^ ^ 

1.036145, and <Pjjj g = .0680697. 

The remaining examples will be derived from the system III.B. 

We shall next consider systems with three environments. Environment 1 
is still the off-environment where no common-cause failure may occur. Environ- 
ments 2 and 3 together constitute the on-environment where common-cause failure 
may occur at a constant rate v*, but with different failure probabilities, 
which we denote respectively by f2 and f^, in order to distinguish from the 
probability q2 in the preceding model. The infinitesimal generator T for envi- 
ronmental changes is chosen as follows : 



T = 



■^12 ^12 

^21 "(^21 * ^23^ 

^21 ^32 



0 

^23 

-(T21 + T32) 



( 6 . 1 ) 



thus the on-environment has exponential duration with parameter T2J, just as 
for the systems III. A and III.B. We shall set T32 > T23 and f2 < fg, so that 
the on-periods are comprised of intervals of relatively low failure probability 
f2, interrupted by shorter intervals of higher failure probability f^. 
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In order to make meaningful comparisons with the systems II. A and 

III.B, we have set, as before, v* = v (1 + expected 

time between successive common-cause failures is the same for all systems. 

We also choose f2 and f^ such that 

.^2 - ^ f 2 ( 1 “^) ^ 3 » 

where a — (^2^ ^ ^ ^32 ^ ^23^ ’ 

and is the conditional probability of being in environment 2, given that the 
system is either in environment 2 or 3. This choice of f2 and f^ ensures that, 

on the average, the probability that any machine fails when conmon-cause fai- 

lure occurs is the same for all systems. 

The system III.C is defined as follows : N = 10, R = 1, K = 3, = X2 

Xj = .005, = 1» 1^2 ~ ^3 * ^2 * '’s * *27666, = 0, f2 = .1504, 

f^ = .36, Tj ^2 “ .00015, T2 j^ = .004, T 22 = *04, T^2 = .125. For system III.D, 
the differences are as follows : f2 = .1801, f^ = .7, T^2 = 1* Thus, we have 
for system III.D very short bursts of high failure probability. 

The stationary probability distributions are respectively given in 
Tables V and VI : in column 0 the marginal distribution of the number of ope- 
rating machines; in columns 1 to 3 the conditional distributions, given the 
environment state. The last row of each table gives the marginal distribution 
of the environment state. 

We compare Tables V and VI with the columns III.B of Table III. The 
same observation that were made earlier are valid here : the marginal distribu- 
tion has a greater probability mass at the states with low operating machines. 
The probabilities corresponding to 6 to 10 operating machines given the off- 
environment (columns III.B.l, III.C.l, III.D.l) are nearly equal in all three 
systems. The conditional distributions given the on-environment are noticeably 
different for all systems; in particular, the distribution in column III.D. 3 is 
different from all the others, in that it is more or less uniform for n between 
2 and 8. 
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In Table VII are given, for the same two systems, the values of 
^10 1 j ^ 3, 0 < k < 9. The last row indicates the conditional 

probability of being in each environment, given that there are 10 machines 
in operation. We observe that short bursts of higher failure probability 
are more detrimental to the performances of the system, since the values of 
UiO |^(j) for the system III.C are nearly consistantly greater than those for 
the system III.D, the exception being u^p j^(3) for high values of k. That 
exception is probably due to the very short duration of the environment 3. 

Finally, we have computed that ipjjj = .0679624, and (pjjj ^ = .0680954. 

Our general conclusion is as follows : if one considers the systems 
III. A (or III.B), II. A and I. A as increasingly simplified approximations of the 
system III.C (or the system III.D), we obtain increasingly optimistic measures 
of their performance, based on the expected first passage times u^q although 
the systems have nearly equal stationary expected number of failures per unit 
of time. 

Type IV . For our last numerical examples, we shall assume that the server does 
not change instantaneously the service rate when there is a change of environ- 
ment. More precisely, we assume that the environment evolves as in the system 
III.C and III.D, with environmental changes governed by a matrix T of the form 
(6.1). The "normal” service rate in the off-environment (environment 1) is 
equal to uj. When the system enters the on-environment (environment states 2 
and 3), the server notices it when 2 or more machines fail simultaneously, or 
after a random interval of time of exponentially distributed duration, with 
parameter Cj. The server then changes the service rate to U 2 - service 
rate remains equal to ^2 until an interval of time has elapsed without any 
failure, that interval of time is exponentially distributed, with parameter ^ 2 - 

In this case, we need to describe the state of the system by three 
indices (n,j,k), where n denotes the number of machines that are operating, j 
denotes the environment state, and k = 1 or 2 to indicate the service rate that 
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is in effect. For each given value of n, we order the states as follows : 
(n,l,l), (n,2,l), (n,3,l), (n,l,2), (n,2,2), (n,3,2). The possible tran- 
sitions from (n,j,k) to (n,j',k') are indicated below, together with the 
corresponding transition rates. 




The infinitesimal generator Q has the form (1.1), with blocks of 
order 6 given by 

P ^ • • • • • 

• ^ X * * * * 

• • ^ X * * * 

• • • ^2 * * 

• • • • ^2 * 

• • • • • p2 



^n,n+l 



, 0<n<N-l, 
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^ 2.1 * 



1.2 • 

T 



2.3 



T T ♦ 
' 2,1 * 3,2 



1.2 



' 2,1 

^ 2,1 ^ 3,2 



' 2,3 



, 0 < n < N, 



nXj 



Q 



n,n-l 





nXj . 





1 < n < N, 



''n.n-k 




q 



( 3 ) 

n,k 





2 < n < N, 
2 < k < n. 



where ,(’j = e v,(J) (1-f,)"-^ 
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the entries represented by a are equal to zero, the entries represented 
by a are such that the row-sums of Q are equal to zero. 

We define the following type IV system : N = 10, R = 1, K = 3, 

* X2 * X^ “ .005, “ 1, ^2 — 2, “ 0, \>2 — ® — .27666, q^ — 0, 

f2 = .1801, f3 = .7, Tj 2 = .00015, = .004, T23 = .04, T32 = 1, ^2 = 

v*/2 = .13833. 

The column 0 of Table VIII gives the marginal distribution of the 
number of operating machines, in the remaining columns are given the condi- 
tional distributions, given (j,k) where j is the environment state, and k the 
index of the service rate. The last row gives the marginal stationary proba- 
bilities of the pairs (j,k). In Table IX are given the values of Ujg j,(j,k), 
for 1 < j ^ 3, k = 1,2 , 0 < n < 9. The last row indicates the conditional 
probabilities of being in the state (10,j,k), given that there are 10 operating 
machines. 

We observe from the last row of Table VIII that the states in {(n,l,2), 

0 < n < 10} are essentially states of transitions towards the set {(n,l,l), 

0 <n < 10). For both j = 2 and 3, the sets {(n,j,l), 0 < n < 10} and 
{(n,j,2), 0 < n < 10} have similar probabilities, but different distributions : 
the states in {(n,j,l), 0 <n < 10} are more concentrated than those in {(n,j,2), 
0 < n < 10} towards large values of n, in other words, towards few machines 
under repair. This is seemingly in contradiction with the inequality < P2; 
this is due to the fact that transitions from (j,l) to (j,2) occur 
whenever there are multiple failures. 

The expected first passage times are not much different from those of 
the system III.D (see Table VII). 



35 . 



Finally, we mention that the program implementing Algorithm E 
has been written in FORTRAN 77 and run on a COC 170/750. Care has been 
taken not to waste memory space, by making sure that no space is reserved 
for matrices n. . with j > i, Q,- .= with j > n+1, or vectors u. . with j > i. 
In other respects, the implementation was straight forward. For systems with 
N = 10, the program uses on the average .2 seconds of central processor time 
for K = 3, and .6 seconds for K = 6. 
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Table I 

Stationary distribution for systems I. A and II. A. 
Missing numbers are less than 10'^. 
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Table II 

Expected time, starting from 10, until there remain 
k or less operating machines, for systems I. A and 
II .A. 
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Table III 



Stationary distributions for systems III. A and III.B. 
Missing numbers are less than 10”^. 
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Expected time, starting from (I0,j), until there 
remain k or less operating machines, for systems 



III. A and III.B 
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Table V 

Stationary distributions for system III.C. 

- 6 

Missing numbers are less than 10” . 
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Table VI 



Stationary distributions for system III.D. 
Missing numbers are less than 10 
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Table VII 



Expected time, starting from (10, j), until there remain 
k or less operating machines, for systems III.C and III.D. 
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Table VIII 



Stationary distribution for a Type IV system. 
Hissing numbers are less than lO'^. 
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Table IX 



Expected time, starting from (10,i,k), until there 
remain n or less operating machines, for a type IV 
system. 
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APPENDIX A 

First order derivatives of the matrices D^(c) 

The equations (3.13) and (3.14) may be written as 

0 q {^){ el - Qq^q) = I, 

- On.n - Dn.k * >• 

By differentiating both sides of these equations with respect to we 
obtain that 



DQ(e)( el “ Qq^o^ ^ ^0^^^ ~ 



n-l 


n-l 


- <*n.n - <!„.k ^ 


<■ - ■ “• 




1 < n < N-l, 


or DJ(5) » - (D„(t))^ , 


(A.l) 


n-l 




o;(e) = - D„(0(I - Q„.k <^k.n(^)) V^)> 


for 1 < n < N-l. (A. 2) 


First order derivatives of the matrices Dp(e) 




The equations (4.5) and (4.6) may be written as 




- Qn,n) = 





- Qn.n - Qn.n.l Vl.n(^)) = Kn<N-l. 

By differentiating both sides of these equations with respect to c. we 
obtain 



Sfj(e)(ei ■ ~ 



A. 2 

or B,;(c) . - (DK(5))^ (A. 3) 

and 0'(5)(5l - Q„.„ - * D„(E)(I - G^„.„(E)) - 0. 

t’(E) = - 0„(5)(I - 5;„i,„(t)) 0„(5). for 1 <n <N-1. (A.4) 



or 



APPENDIX B 



Algorithm C . Evaluation of the stationary probability vector £ = {£o» 
and the expected first passage times u_ , for all pairs of levels n 

Cl (Initial computations : the matrices |^, and the vectors ^). 

^N ’ 

determine Uj^ from (4.12); 

for k from 0 ^ N-1 ^ determi ne rij^ from (2.4) ; 

for n from N-1 1 step -1 do 

begin determine from (2.3); 
determine u^ from (4.13); 

for k from 0 n-1 ^ determi ne ^ from (2.5) 

end; 

determi ne Cg from (2.3) ; 

C2 (Determination of the vectors £^ and u^ ^,) 

iig := solution of the system £q Cq = 0^» £o £ = li 

^0 %,0‘ 

- 1,0 •“ - 1 ’ 

for n from 1 ^ N ^ determine ^ g from (4.16) ; 
determine Ug ^ from (3.17); 
determi ne Fg j from (3.4); 
for m from 1 to N-1 do 



B.2 



begi n determine from (2.7); 

m 

re-normal ize the vectors to ir„ so that Z w. e = 1; 
-0 i=o - 

for n from nH-1 to N do determine u from (4.16); 

— n,m ' 

determi ne from (3.3) j 

determine from (3.18); 

determine from (3.4); 

for n from 0 ^ m-1 do 
begin determine from (3.19); 

determine from (3.5) 

end 

end ; 

determine Hm from (2.7); 

^ N 

renormalize the vectors to tij, so that Z n„ e = 1. 
-0 -N n=0 - 



At the end of the algorithm, the vectors are equal to the vectors £^, 
0 < n < N. We systematically renormalize the vectors in order to avoid 
overflow problems. 



Algorithm D. Recursive evaluation of the vectors defined ^ ~ » 

0 < k < N-1. 

D1 (Initial computations) 

determijTe ^ from (4.12); 

^-1 (equation (4.17)) 



B.3 



02 (Main loop, at the k th step, is determined) 
for k from N-1 ^ 1 step -1 ^ 

begin for n from N to k+1 step -1 do determine n_ from (2.4) and (2.5); 
determine C|^ from (2.3); 
determine ^ from (4.13); 

^,k-l- '" (equation (4.17)) 

for n from k+1 N ^ determi ne ^ |^_j^ from (4.18) ; 

^-1 H^.k-i 

end . 

Algorithm E . Evaluation of the stationary probability vector £ = (pg g^), 

and the expected first passage times ^ for 0 < k < n < N. This algorithm 
combines Algorithm D and part of Algorithm C. 

El (Initial computations) 

^N,N ’ 

determine U|,j from (4.12); 

for k from 0 ^ N-1 ^ determine n.^ from (2.4) ; 

E2 (Determination of the vectors ^ . and preparation to the evaluation of the 

““Tl f K 

vector £) 

^ i (equation (4.17)) 
for n from N-1 1 step -1 do 

begin determine from (2.3); 
determine u^ from (4.13); 

^,n-l i ’ (equation (4.17)) 

for m from n+1 N ^ determine ^ ^ from (4.18) ; 

for k from 0 to n-1 ^ determine from (2.5) 



end; 



E3 (Determination of the vector £) 



tiq := solution of the system JIq Cq * 0.» Hq £ * 1; 

for m from 1 ^ N-1 do 

begin determine from (2.7); 

m 

renormalize the vectors Tin to tt„ so that Z 
-0 -m -1 



end; 

determine Wuj'from (2.7); 
renormalize the vectors iiq 



N 

to tTm so that Z n. e = 1. 

i=o - 



(The vectors are equal to the vectros 0 < n < N) 



C.l 



APPENDIX C 

Infinitesimal generator for Type III . 

We assume N = 4, R = 2 and K = 3. Entries represented by a 
are equal to zero, entries represented by a are equal to minus the sum 
of the off-diagonal entries on the same row. 




where q^*^) = nx^ + (J) qj (l-qk)"’\ 1 < k < 3, 1 < n < 4, 1 < i < n. 
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